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ABSTRACT 

A study is made of the behavior of chemically generated waves in a simpli- 
fied atmosphere. The atmosphere is assumed unbounded, isothermal, one- 
space-dimensional and initially quiescent. At an initial time a dissociation re- 
action, AB + J->A + B + J, commences and drives the subsequent wave motion. 
The fraction of reactant in the atmosphere, X 0 , is assumed to be small. The 
system of governing equations is then expanded in terms of the small parameter, 
X 0 , and an asymptotic integral solution as X 0 -► 0 is obtained. 

The solution is integrated numerically for several natural impulsive heat 
releases, including thunderstorms, hurricanes, chemical reactions, and 
aurorae. 

Strong thunderstorms and hurricanes are found to induce considerable ver- 
tical motion and large temperature changes in the upper atmosphere. The same 
is true for auroral heat releases. Conversely, upper atmosphere effects of 
changes in the ozone layer appear to be small. 
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A LINEARIZED SOLUTION TO IMPULSIVE LATENT HEAT 
RELEASE IN A DILUTE, ISOTHERMAL ATMOSPHERE 

I. INTRODUCTION 

Waves in a stratified fluid under the influence of gravity appear to have 
been initially discussed by Burnside (1889) and Love (1891). Both authors treated 
an incompressible fluid. Gortler (1943) used schlierin photography to show 
experimentally that disturbances in an incompressible stratified medium under 
the influence of gravity propagate along characteristic rays. Lamb (1908) 
treated a compressible, adiabatic, ideal gas whose density is stratified by gravity. 
Since the original contributions by Burnside, Love and Lamb there have been many 
theoretical papers on various types of gravity waves. Reasonably up-to-date 
treatments of the subject may be found in Eckart (1960) and Yih (1965). However, 
gravity wave theory is currently undergoing an active phase of development 
and contributions are being added to the literature at a rapid rate. 

Recent experimental data from the earth's atmosphere and oceans have 
greatly stimulated interest in gravity- acoustic waves. Gossard (1962) has ob- 
served gravity waves in the troposphere. Hines (1960) has shown that internal 
atmospheric gravity waves may account for many of the phenomena observed 
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in tho lower Ionosphere, ami gravity waves have been observed in the thermo- 
sphere by Newton et. al., (1969) and Harris et al., (1969). 

It seems to be generally agreed that gravity waves are generated in the 
troposphere, then propagate upwards. Thus Gossard (1962) has observed 
gravity waves near the earth's surface and Eberstein (1970) has illustrated 
the development of gravity waves between 30 km and 120 km. 

The simpler forms of gravity wave theory treat propagat'on through an 
inviscid, non-conducting gas, while more sophisticated theories include effects 
of heat transfer, viscosity, and high altitude phenomena such as ion drag. How- 
ever, the effect of chemical reactions does not seem to have been adequately 
considered to date. 

Between the troposphere and the thermosphere there are several regions 
where important chemical reactions take place. 

First, there is the ozonosphere between approximately 15 km and 45 km, 
with peak ozone concentration at about 35 km. Ozone mole fraction is in the 
order of parts per million (Mitra, 1952). 

Between 60 km and 96 km there is Nitric Oxide (Pearce, 1969). The mixing 
ratio for Nitric Oxide is also in ppm. 
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Between 90kra and 120 km molecular oxygen dissociates into atomic oxygen. 
The dissociation ratio, a being 3 x io' N at 90km and 0.998 at 120km (Mitra, 
1952). The oxygen is no longer a truly dilute reactant since its mixing ratio is 
20%. However, the high mixing ratio of oxygen becomes a problem only where 
the final fraction dissociated becomes large. Otherwise, one may use the math- 
ematical expedient of considering only those oxygen molecules which ultimately 
will dissociate. 

The major chemical systems are tied to the absorption of ultraviolet radia- 
tion from the sun, and thus have a dawn and dusk dependence. There is also a 
dependence on solar activity. The latter dependence is more important for the 
present study, since we are primarily concerned with impulsive heat release or 
absorption. Periodic thermal effects are in the realm of tidal theory which has 
been extensively treated in the literature. 

Large storms and such phenomena as aurorae may also be simulated by a 
chemical type heat release. 

In what follows , a study is made of the effects on the atmosphere of impul- 
sive heat release in a limited altitude range in a dilute, isothermal atmosphere. 
A chemical reaction is used as the heat release mechanism, and simplifying 
assumptions are introduced to make the mathematics more tractable, and per- 
mit an anlytical solution. Our objective is a relatively simple theory which may 
be readily used to study the overall atmosphere effects of large classes of 
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natural impulsive heat releases. It is recognized that comprehensive numeri- 
cal techniques have been developed to study the atmospheric effects of impulsive 
heat release from a point source. However, complicated computer programs 
are neither generally accessible nor readily modified for general scientific use. 
Thus we feel that a simple theory of the effects of impulsive heat release in the 
atmosphere could be helpful in spite of limitations imposed by simplifying 
assumptions. 


II. GOVERNING EQUATIONS 

The continuity and momentum equations (cf. Shere and Bowhill, 1969) are 
not affected by the chemical reaction provided that the state variables are given 
the extended definition: 


I . P(p, T. a) 


wnere P, p t T and a denote pressure, density, temperature and degree of 
advancement of the reaction, respectively. These equations are in tensor 
notation: 


D. * u i 

m + p TT = 0 

i 


( 2 . 1 ) 


and 


Du, ,?p da. t t _ 

p Dt ' n , * ax, + a Xj (* = 2 > (2.2) 


where 


a . 

) » 



7*7 


<* u j 2 du k~ 

7Tx~ T i j e)x k 
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We denote time, space variable (x ; vertical), velocity component, gravitation 
vector component, viscous stress tensor, dynamic viscosity and the Kronecker 
delta by t, x ( , u ( , g ( ( , m and respectively. We also take g ( «-g& |2 
and we define D/Dt to be „he total or Eulerian derivative. 


In order to develop the energy equation, it is necessary to consider the 
thermodynamics of the system. The derivation of the equation of state parallels 
the work of Liepmann and Roshko (1957, p. 29) and the derivation of the reaction 
rate equation parallels the work of Eber stein (1966). The following discussion 
applies to the reaction 


J + AB s: A + B + J 

where J is an inert element. We assume that each component of the mixture is 
a perfect gas. The mixture, however, is not a perfect gas. The partial pressure 
of the ith constituent is 


P, 


m 


. P R, T 


(2.3) 


where R f K/W ( with K the universal gas constant, W t the molecular weight, 
and m ( the mass fraction of ith-component. Letting a be the fraction of AB 
dissociated, n ( the number of moles of ith component, and n 0 the total number 
of moles when a = 0, 


n A B 

o 

ll 

x 0 (1 a) 

n A 

o 

C 

II 

X o a 

J 

00 

ll 

mJ 

o 

X o a 

n J 

r (1 

X 0 ) n o 
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where X 0 is the mole fraction of AB when » 0. Summing yields 

n. n 0 ( 1 ♦ a X 0 ). 

P = K( a)//T (2.5) 

R( a) = i M 


where R is the universal gas constant, and M is the mean molecular weight, 
given as: 

M = [*j + X 0 (W A „ - *)]/[l ♦ a X„] 

It follow s that 


R(a) = R( 0) (1 + a X 0 ) 

The reaction equation is as follows: 


D "AB /Dt 



(2.6) 


(2.7) 


where 



and k is a rate constant. In determining (2.7) it has been implicity assumed that 
all molecules present may act as an ineit element with equal collision effective- 
ness. Substituting (2.4) into (2.7) yields 


Dt 


[pk F fT) M] 



U 


[K(T)/o/M] 



( 2 . 8 ) 


where 


Igor J. E be rate In 


and M is the mean molecular weight 


For chemical equilibrium, the rate equation reduces to 



where use was made of the equality 



In atmospheric problems the reactant fraction is generally small, and it 
will be assumed that X n << 1. It follows that aX Q << 1. The rate equation 
may thus be simplified to give: 

dT = [^f<T)/m] {l - a} (2.9 

The energy equation is now developed in the usual manner (cf. Liepmann 
and Roshko, 1957, p. 185+). Summing the rate of change of the internal, kinetic 
and potential energy yields 
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where n is the ith component of the unit normal of the closed surface />, V is 
the volume of i and A the area; e is the internal energy per unit mass and * is 
the coefficient of thermal conductivity. 

By applying Gauss’ theorem and substituting the continuity and momentum 
equations into (2.10), we obtain 



We now use the first law' of thermodynamics, h - e + P/p, and the continuity 
equation to obtain 



In a reacting gas mixture the enthalpy depends on the temperature and the 
degree of reaction, i.e. 

h = h (T, a) 

The rate of change of enthalpy per unit mass is then given by: 


Dli DT Da 

■ET = c ,¥ ,x o B Dt 


( 2 . 13 ) 


where B is the enthalpy of reaction and C p is heat capacity. For an endothermic 
(exothermic) reaction, B is positive (negative). Both the enthalpy of reaction 
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and the heat capacity have a weak temperature dependence, and an even weaker 
pressure dependence. For simplicity both C () and B will be assumed constant 
throughout this paper. 

Eliminating h between (2.12) and (2.13) one obtains the enthalpy equation 
in the form: 



♦ X 0 B 


C a 
Dt 



(2.14) 


The independent variables in the system of equations described above 
are t, x . The quantities k K , K, C ( and B are known parameters of the system, 

X 0 is a small (known) parameter. The dependent variables are P, p, T, u , and a. 

For the remainder of this paper we assume that the atmosphere is inviscid 
and nonconducting, i.e. , * = X =0 

HI. NONDIME N SIONA LIZ AT ION OF THE GOVERNING EQUATIONS 

Since mathematical operations are performed on pure numbers it is proper 
to non-dimensionalize the governing equations. 

Some of the non-dimensionalization is almost trivial, as shown immediately 
below: 


R' 

- R R(0) 

M' 

= M/M* 

C ' 

= C p R(0) 

V 

k F P m 

9 

P 

= p/p . 

K' 

= K M*o; b 
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where H(0) is the gas constant for > = 0, and the starred reference level is the 
value at the bottom of the reaction zone. < B is a base frequency, and H is 
scale height 


t' 

- t 

1 

X . 

= X >1 


B 

1 

1 

T' 

= R(0)T/ &B 2 H J 

I 

g 

r R H' 

B' 

= ■V * 1 

P' 

p f p 

* 


It is ? 0 en that in the primed system we have 

P' p' R' T # 

In general, the non-dimensional system of governing equations looks very 
similar to the dimensional system. 


The quantity < B 2 H 2 is proportional to the energy contained in a wave 
oscillating with the base frequency, and having amplitude equal to a pressure 
scale height. The quantity H B 2 is a measure of the acceleration experienced 
in such a wave, ^ince such acceleration must be less than that due to gravity 
we would expect g' > 1. 


One may estimate the magnitudes of the non-dimensional parameters. 
Taking a scale height of 6 km and a base period of 5 minutes, one obtains 

(a 0 H) 2 ~ 1.6 * 1 0 8 cm 2 /soc 2 3.74 calories/gram 

and 

_ _ # g ' ~ 3. 5 

cV H ^ 261 cm/sec 2 , giving 

T ^ t 4 ( 
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Further, taking B = 34.4 kcal/gmole, M = 48g/gmole, one obtains B* ' 200. 

IV. ASYMPTOTIC DEVELOPMENT 

The non-dimensionalized system of equations will be expanded about the 
parameter X 0 , i.e., we expand each dependent variable in a power series of 
the form: 


f ( t , x , z ) 


L 

N - 0 


f c N ) ( t . X, z ) 


X N 
A o 


(4.1) 


where higher order terms may be dropped as X 0 *0. r i'his expansion limits 
the validity of our results to dilute fluids, such as the atmosphere. Other ex- 
pansions [cf. Lax, 19571 can be obtained by inserting a small parameter into 
the initial conditions. These expansions may, for example, be used to study 
the generation of gravity waves by tidal waves. 


Assume that the atmosphere is initially stratified, quiescent, and bounded 
below by the earth, i.e., we define an initial-bourn i.y value problem with initial 
conditions: 


T( 0, x , z ) 

1^(0, x, z ) 

P( 0, x , z ) 
a( 0, x , z ) 


T - T 
o , 


0 


e 


(i = 1, 2) 


0 


(4.2) 


and boundary condition: 
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u 2 ( t , X , 0) = 0 

Using the equation of state to eliminate P from the governing equations and 
substituting (4. 1) into the resultant system yields the zero-order solution 


u 


T< o > 

(t. 

x, z ) 

T = T 
1 o 1 . 

(0) 

i 

(t. 

x, z ) = 

0 

p<0> 

(t, 

X, z ) = 

e~ * 

a<°> 

(t. 

X, z ) = 

1 ~ exp [**k F 0 * t 

rj. ( 0 ) 

T • p( 1 > 
A d) * 

= p ( 0 ) p^ X) and a ' 1 > 


(4.3) 


(I ) 1 


we get the first order system 


>u<‘> 

dt 


+ T 


' ( i ) 

o x 




» T 


( i ) 


0 (1 x 


0 


(4.4) 


>u a < 1 > 
fit 


it . 


< i ) 


0 r)z 


^^(i) x Ba (0) 

4 T — - — L - T, 


(0 ) 1 ( 1 ) 0 dz 


0 3z 

+ T 0 a <0) 


< 1 ) (3m U ) 

( l ) U l 0U 2 

— wr \ k 4 3 u„ ( 1 ’ 0 

at i7X a z 2 


(4.5) 


(4.6) 


d P<n 1 dT d ) 

d t * y ~ 1 <5 1 


+ 


3 aa(°) 

~ 9 IT” 


(4.7) 
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(°> ♦ a 


,1 i<° > 


( i ) <9t 


= T <n (p<°>, a<°> ) 


(4.8) 


where 


8 = 


= (t» ' ')• 


Equations (4.4) - (4.7) form a linear system of four partial differential 
equations in four unknown functions. The system is inhomogeneous. The solu- 
tion of (4.8) for * is needed only for the calculation of second order terms. 
Since we determine only the first order terms in this paper, we will neither 
specify nor solve for o. (J ( . 


Dropping the (1) subscript or superscript and eliminating we obtain the 


d 2 [dj 

dt 2 [at 


2 * 2 7 d 2 

- c 2 — T “ c 2 — T 

dx 2 d z 2 


r] h «• 


where 


h - (r- 1 ) 


8 it T L^t 2 t, j, 


2 d 2 

_ - T — a<°> 

.2 0 •, „ 2 


‘ (y_ 1 )B T 7 ^T + (7 “ 1)T » 
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and 


C 2 



After (4.9) is solved, the other dependent variables can be determined by 
the equations 


r "■ . t ■ ] 

n,/ J 


— . - T "T - + T T + 8~ — ~ + T. a<°> (4.H) 

- 1 . * n * Vt 2 ° 


I_ d 2 T 81 

7 “ 1 ^ t 2 ^0 d Z *0 


dll 


s 2 . JL '* 2 t JL ( , Jifi. 

d t d z + 7 - 1 d t 3z 7 “ 1 d t dt\dz 


- l) a<° 


(4.12) 


and 


u. 


dp 1 dT da<°> 
TT " T^T 7TT “ B 


(4.13) 


Now define 


#(t, x, z) * e~ ‘ 2 T(t, x, z) 


(4.14) 


The above definition allows the temperature perturbation to grow expo- 
nexitially with altitude without a corresponding growth in 6 . 

Substituting the definition of Q into (4.9) yields: 


_d^ 

dt 2 



C 2 V 2 


c 2 

+ -4- 




e - 


C 2 (7“ 1) d 2 6 


r 


dx' 


e 


-z /2 


f (a) 


(4.15) 


f(a (0) ) was defined in Equation (4.10). 
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The dispersion equation of (4.15) is: 

V - c * (\ 2 + x , 2 , t) X , 2 * [(y- 1 ) C*/y] x, j = o ( 4 . 16 ) 

where Is the wave number of i. 

For a discussion of dispersion relations and how these are obtained from 
partial differential equations reference is made to Courant anti Hilbert, es- 
pecially p. 588. 

V. SOLUTION OF THE FIRST-ORDER SYSTEM 

From (4.4)-(4.6) and the initial conditions (4.2) we deduce that 


0(0, X, 

z) - 0 

(5.1a) 

0 t (O. x, 

, z) = - (y - 1) 8k f e" 3,/2 

(5.1b) 


Conditions for 0 tl and & %xx at t = 0 can be determined by differentiating (4.4)- 
(4.8) appropriately and analyzing the resultant system or by integrating (4.15) 
and applying the physical requirement that £(t, x, z) must remain bounded as 
t -oo. The negativeness of 5> t (0, x, z) means that the atmosphere initially 
cools, which is in correspondence with the endothermic reaction indicated by 
positive 8 . 

Since the inhomogeneous parts of (4. 5) -(5. 1) do not depend upon the hori- 
zontal space variable, neither does the solution. Thus (4.15) reduces to 
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Mf/?1 


* / 2 


f( o 


(5.2) 


where 


Mi*) d 2 /<? 1 2 ”C 2 d 2 /<? i. 2 + C 2 /4 . 


(5.3) 


It is worth noting that the initial conditions and the alpha dependence chosen 
are such as to imply the relevance of only one space dimension. Other initial 
or boundary conditions would require that the second space dimension be kept. 
Unfortunately, exclusion of the lateral space dimension removes many interest- 
ing physical phenomena from the realm of applicability of the theory. However, 
it must also lx? realized that keeping the lateral space dimension is only mean- 
ingful if one also knows the horizontal behavior of the forcing function. Other- 
wise, one can only calculate propagation of waves, but not their generation. 


Integrating (5.2) with respect to time twice, 

M[tf] - w ( t , z) = k F (y - 1) ik f [-8 + Bt 2 ] e"* - 2Btlcxp 
subject to initial conditions (5. 1). 


3z 

T 


+ kpte 


■] 


(5.4) 


The inhomogeneous part of (5.2) is the effect of the chemical reaction. For 
a realistic approximation to the atmosphere, the inhomogeneous part must 
be multiplied by a suitable weighting function, since the chemical contribution 
only extends over a finite altitude regime. We also note that at each fixed 
altitude the driving force must tend asymptotically to zero in the Poincare 
sense as t - , . 
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The solution of (5.1)-(5.4) is obtained by superposition. We choose a suitable 
function ] (t, z) which satisfies the initial conditions and is damped in time. 

For example, set 

' j (t , z) - (y - 1 ) ' kj t e” 2 ,c e~ 3 * 2 (5.5) 

Observe that (5.1) is satisfied and 

MtOj] = 2/2 C (y - 1) Bk f e‘‘ 2 tc e*“ J,/2 

Since Ml ] is a linear operator, (5.1)-(5.4) may be reformulr.ted as 

0(t, z) = ^(t, z) f 0 3 (t , z) (5.6) 

where < • (t, z) satisfies 

M[0 a (t, *)] -w(t. z) - 2/2 C(y - 1) 8k f e~ ltc e’ 3l/2 f (5.7a) 

and 

0 2 (0 . z) r e 2 t (0, z) = 0. (5.7b) 

The solution of (5.7) for t * 2 (t, z) is known Icf. Lamb, 19091 and is given by 

& 2 (t, z) = f W(t, z; r) dr (5.8) 

•'o 

where 


W ( t . z; r) 



• Q(r, z, t?) dr? (5.9) 


17 


Igor J. Eberstein 


and 


with 


Q(7, *, Tj) [ i (r, z * r ( ) ♦ >p(r, z - rj)] /2c 2 


g ( • , z ) w ( - , z ) - 2 7c(y - 1 ) - k f c~ 


(5.10) 


Although (5.5) and (5.8)-(5.10) provide an exact solution of (5.1)-(5.4), this solu- 
tion is complicated and yields quantitative information only through numerical 
calculations. The results of these calculations are given in the following sections. 


A solution of (5.1) -(5.4) can also bo obtained using series techniques. The 
result is a somewhat rapidly convergent infinite series whose terms Involve 
products of powers (- k ( te ')" and Young’s functions Young, 1912 , a special 
case of Lommel functions. Since the series was not used to obtain physical 
results, it is not presenteu here. 

From (4.12)-(4.13) we deduce the first order terms: 

T( t , z ) = c* 2 0 ( t , z ) 


u 2 (t. z) 



+ 



dz 


where z is at the earth surface and 
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VI. COMPUTING 

The cjre of the numerical problem is evaluation of the integral expressions 
(5. 8) and (5. 9) 

Evaluation of the above integrals was done numerically using Euler’s 
method : 


N 

£(t,z)= y W t (t, i. r. ) At 

I . 0 


(«. 1) 


where At t/N and r. iAr 
Similarly, 


W 


yi 

■ L 


i • 0 i— 


J o ll - T i > 2 - (t 


Q( t |. z. V } ) At7 ( 6 «-) 


where 

At;' (t - r) c/M 
*?; = j At; 

The Bessel function was evaluated by numerical integration to at least five 
significant figures. Accuracy o f the routine was checked against tables in the 
handbook of Chemistry and Physics. 

The size of the integration mesh was decreased until the results became 
invariant to further decreases of mesh size. 
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The above procedure works quite well for short or medium times. If the 
integration steps are kept small for longer times, then computing becomes 
excessive. Conversely, if the number of intervals is kept constant, then 
accuracy suffers. Also, it is desirable to observe how the perturbation profile 
develops with time. We thus break up the integral for 0 2 into a series of inte- 
grals as follows: 

>L, /«Atk 

0 2 (t,z)= > W(t k ,z.r)dr (6.3) 

kt| •'At(k-I) 

where K is chosen such that 

t = KAt 


Since W is a function of final time, (6.3) cannot be solved unless an addi- 
tional boundary condition is imposed. 

Let us examine (5.9): 

The quantity \l (tc - tc ) 2 - rj 2 is familiar from the theory of wave propaga- 
tion. (t - r)c is the distance that a wavelet has travelled in time (t - r). 

Thus the r, integration goes to c (t - r), or to the limits of the physical region 
being considered, provided that there is no reflection at the boundaries. Let 
emax be the distance from the furthest point in the region considered to the 
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furthest point in the reaction zone. If tc > emax, one has summed all the con- 
tributions. Defining At > emax/c, it is thus permissible to write: 


6L = 


/» At 

* aAt 

r 

W dr + 

W d r + • • • 


Jo • 

Ja» 

J< 


.nAt 


(6.4) 


with 


W 


r 

»'0 


(k At -r ) 


J 0 0 dr? 


(6.5) 


and 

(k - 1 ) A t £ r < k At. 

The quantity t - r ranges from zero to At making the values of the Bessel 
function independent of k. Since the range of rj is thus tc, it becomes possible 
to compute a matrix of Bessel functions 

MJ 0 (I, J) 

where the index I refers to (t - r) (I) and the index J refers to r? (J). Q, however, 
depends on the actual time elapsed, and must thus be computed for eaca value of 
k. 


The integral for W is therefore approximated by: 


J MAX 

W (K, I) = Y M J 0 (I, J) Q (K. J) A 7)(J) 

J * 1 


(6.6) 
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«,=■ V ®i<o 

I « 1 
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(6.7) 


and 


«,(D 


KMAX 


L 


W(K. I) Ar (K) 


(6.8) 


The rate of change of temperature, or 6 % is needed to evaluate the velocity. 
This is computed using the forward difference approximation: 


d G ^ A 6? _ 'n * rn (6.9) 

d t * A t t - t 

n m 


The derivatives of 6 x and a are evaluated from exact analytical expressions. 
The velocity is estimated by a simple numerical integration in z, taking the 
velocity at the earth’s surface to be zero. Density deviation is estimated the 
same way. It should be noted that the step sizes for z were between 0. 1H and 
0.5H. Thus, the velocity and density estimates are generally less reliable than 
the temperature estimates. 

Vn. NUMERICAL RESULTS 

At very short times a pulse is seen propagating up and down from the re- 
action zone. The initial pulse has the appearance of a discontinuity. Eventually 
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the pulse passes outside the range of the computation regime, and a pseudo- 
steady pattern is established in which the qualitative behavior of the parameters 
does not change. However, the quantitative values increase to a maximum, and 
then decay. The development and decay of the pseudo-steady patterns are illus- 
trated in Figures 1 through 3. 

The type of behavior observed may be partially explained by analogy with a 
shock tube whose driven end is semi-infinite. Initially the shock passes, then 
a pseudo-steady state is established and eventually decays. The above analogy 
is incomplete, since the atmosphere also behaves like an elastic medium result- 
ing in the establishment of something like a standing wave pattern. However, an 
acoustic treatment w’ould be incorrect because the gravity restoring force is 
important, and because the reaction generates pulses all of one sign, either 
compression or rarefaction. Perhaps another analogy is a spring with weights 
at the end. When weights are suddenly added or removed, wave patterns are 
set up in the spring. Now, consider that a series of weights are added or re- 
moved in succession. Also, let the spring be very stiff near the bottom, becom- 
ing progressively more elastic toward the top. 

Figure 1 shows a series of 0 profiles. The length of the reaction zone is 
0.5H, the reaction is exothermic, and the nondimens ional rate constant is 0.01. 
Figure 1 shows buildup of the 6 profile. The reaction zone was chosen to be 
0.5 H deep, the reaction was exothermic, and the nondimensional rate constant 


23 


Igor J. Eberstein 


w'sO.Ol. The reference value, K , equals lO^sec -1 ; thus time is given in units 
of 100 seconds. The characteristic time of the reaction is defined as the point 
where k F t = 1. Strictly speaking, one needs k F ^ (0) t - 1, but /o (0 > = 1 at the 
bottom of the reaction zone. Since ; (0) drops off exponentially with altitude, 
while k F remains constant, it follows that the characteristic time increases ex- 
ponentially with altitude until the end of the reaction zone is reached. It may be 
shown that 


k F t = k F * t* 

where the starred quantities are dimensional. For k F = 10" 2 sec -1 the charac- 
teristic time is 10 2 in non-dimensional units, or 10 4 sec = 2.8 hrs. At z = 0.5 H, 
the characteristic time becomes 4.G hrs. 

The absolute maximum value of 0 occurs at z = 0 and grows with time. A 
relative maximum is found at z = 10, this secondary maximum al r o grows with 
time. Essentially, 6 follows a Bessel function type of altitude pattern as might 
be expected. The growth of 0 at z = 0 and at z = 10 is presented in Figures 
2 and 3. Figure 4 shows development of the temperature profiles. 

Outside the reaction zone the velocity equation becomes: 

Bu ~ 1 3T 

B z y - 1 3 1 


During the buildup of temperature we have 
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everywhere except in a small region near six scale heights. One may conclude 
that velocity is generally negative during the buildup phase, especially above 
8 scale heights. This conclusion is indeed correct. 

Physically we know that subsidence in an isothermal atmosphere results in 
heating. Conversely, an upward motion gives rise to cooling. This is precisely 
the kind of behavior observed, so we may conclude that the result is physically 
consistent. Similarly, convergency of velocity results in compression while 
divergence gives rarefaction. These combined effects are illustrated in Figure 
5 which, incidentally is for an endothermic reaction, five scale heights deep. 

If we had an adiabatic lapse rate, then physical argument would lead us to 
expect disappearance of the most pronounced part of the temperature wave. 
However, the velocity wave would not disappear; consequently, a temperature 
effect would become visible at an altitude where the lapse rate became less than 
adiabatic. If the atmosphere has an inversion, then the temperature wave would 
be especially pronounced. 

Several special cases are discussed below. 

The first case to be considered is a severe rain storm. Condense saturated 
air at 5km altitude (270 degrees Kelvin) imtil 25% of it is ice, using a reaction 
characteristic time of 3 hours and a reaction depth of half of a scale height. 

The mixing ratio, X Q is 0.002, and the non-dimensional enthalpy, B, is just un- 
der 200. The reaction is of course exothermic. Since the extent of a rainstorm, 
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or even a hurricane, is not large enough to really warrant a one-dimensional 
approximation, the storm was assumed to be approximated by a line source 
centered at G II below ground level. The above approximation results in a ■— 
amplitude fall-off relative to the pure one dimensional case. It is realized that 
the two-dimensionality correction employed is quite arbitrary, and thatagenuine 
two-dimensional solution of the atmospheric equations of motion is needed. 
However, the results, shown in Figure 4, seem to agree reasonably well with 
experimental observations, as shown by Eberstein and Theon (1971). 

Our calculations were taken to ten scale heights above the bottom of the 
reaction zone. Whereas the ten scale height limit was chosen quite arbitrarily, 
there are nonetheless compelling reasons for limiting the vertical extent to 
which computations are carried out. Firstly, the one-dimensional assumption 
becomes ever less meaningful as the vertical extent of space is increased. 
Secondly, the temperature perturbation involves an exponential in altitude, i.e. , 

1-6 e ‘ /2 

with the consequence that small errors in 6 can give rise to large temperature 
errors as z becomes large. Also, the non-dissipation and isothermal atmos- 
phere assumptions lose validity as one considers effects propagating over large 
distances. A more detailed and comprehensive theory is needed to study the 
effect of severe thunderstorms or hurricanes on regions in the ionosphere and 
above. However, the present Jinalyses illustrates that such effects would 
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definitely Ijc expected. This conclusion is borne out by experimental observa- 
tions. Thus Bauer (1958) has shown a correlation between hurricane passage 
and electron concentration in the F 2 layer of the ionosphere. More recently, 
Davies and Jones (1971) have reported association between ionospheric distur- 
bances in the F, region and severe thunderstorms. Davies and Jones believe 
that the ionosphere is perturbed by infrasonic disturbances generated by mechani 
cal motions of the thunderstorms, but are not due to buoyancy oscillations. We 
would suggest that heat released by the storm induces vertical motion, thus 
influencing the electron concentration and transmission properties of the F 2 
layer. 

While the ionosphere is outside our quantative reach, the ozonosphere is 
relatively accessible. Reed (1950) suggests qualitative explanations in terms of 
vertical and horizontal motions for correlations between ozone concentration 
and weather phenomena. The existence of an ozone-weather relationship is 
described as well known. Reed specifically considers subsidence at high alti- 
tudes as one of the means by which ozone concentration is increased. Our model 
predicts that a severe storm will cause considerable subsidence at ozone alti- 
tudes, in agreement with Reed’s suggestion. 

At the time of computation the velocity and density information were con- 
sidered to be of secondary importance. A rather crude z mesh was thus used 
to save computer time, with the result that the quantative velocity profiles must 
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be considered approximate. Nonetheless, we are confident that a severe storm 
causes large and sustained subsidence at ozone altitudes. It would be very in- 
teresting to use a fine vertical mesh computation or a more sophisticated mathe- 
matical technique to obtain a quantitative estimate of the ozone concentration 
change. 

One might also consider the effect of changes in the ozone layer on the 
upper atmosphere. A reaction having the thermal properties of ozone dissocia- 
ting to molecular oxygen was considered. The reaction characteristic time was 
taken to be 20 minutes. This results in a maximum perturbation at approxi- 
mately 4 hours. The type of perturbation profile attained is shown in Figure 5. 

A 2 degree Kelvin cooling at the first maximum seems quite reasonable (Krueger, 
1971). At the second maximum (911 above the bottom of the reaction zone) the 
temperature change is 15 degrees. The associated density perturbation is 6%, 
and vertical velocity is 40cm/sec at 10 minutes, going to 5 cm/sec at 4 hours. 
Since the second maximum is above 90 km where large atmosphere variations 
are frequently found, one must conclude that upper atmosphere effects of ozone 
variations are not very important. 

A very different conclusion is reached regarding the effects of aurorae on 
the upper atmosphere. Insofar as a rather large amount of heat is rapidly re- 
leased in a small altitude regime, the thermodynamic effect of aurorae is very 
similar to that of a severe thunderstorm or hurricane. The expected behavior 
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is thus generally similar to that shown in Figure 4. If one takes an initial heat- 
ing rate of 25 ergs /cm ‘ -sec and a reaction characteristic time of 20 minutes, 
with heating concentrated in half a scale height, then the quantitative deviations 
are approximately the same as those snown in Figure 4. If the mean heating 
rate is 25 ergs/cm 2 -sec, then the deviations become twice as large. According 
to D. Heath (1971) auroral heating rates vary between 10 ergs/cm 2 -sec and 100 
ergs/cm 2 -sec with heating concentrated in less than one scale height. Corre- 
sponding characteristic times vary between 10 minuteR and 100 minutes. Since 
our solution is linear, it becomes possible to estimate upper atmosphere effects 
anywhere in this range. The maximum temperature deviation predicted is then 
in the order of 500 degrees at some 200 km for the case of a mean heating rate 
of lOOergs/cm 2 -sec. The actual value of 500 degrees must, of course, not be 
taken too seriously, especially since large heat releases, and associated large 
accelerations violate some of the assumptions on which our solution is based. 
Nonetheless, the theory does predict a large temperature increase well above 
the main auroral display altitude. 

vm. CONCLUSIONS 

A one-dimensional model for impulsive heat release in the atmosphere has 
been developed. The theory described is intended as a simple tool to study the 
effects of impulsive heat release. Such heating, or cooling, is found to cause 
large disturbances at higher altitudes. 


29 


Igor J. Kberstein 


The inclusion of a second space variable and additional mathematical model- 
ing of heat sources is needed. Among existing solutions including more than one 
space dimension are acoustic waves, internal gravity waves, and tidal w r aves. 
However, all three above mentioned waves are social restricted solutions to 
the atmospheric equations. Specifically , acoustic and gravPy waves have sinus- 
oidal space and time behavior. Tidal waves have a sinusoidal time behavior 
and a spatial behavior described in terms of Hough functions. The above theories 

are quite good for evaluating the long distance propagation of periodic distur- 
bances. However, these theories may not readily be employed to study the short 
distance effects of impulsive heat releases in the atmosphere. Our one- 
dimensional theory has been an initial step toward an analytical solution to the 
problem of impulsive heat releases in the atmosphere. Thunderstorms, hurri- 
canes. chemical reactions, and aurorae have been discussed ? i important natural 
sources of impulsive heat release. 


' 
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TIME IN 100 SECOND UNITS 


Figure 3. #* Development with Time at z - 0 
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ALTITUDE; SCALE HEIGHTS 


1500 SEC 



TEMPERATURE DEVIATION; DEG. KELVIN 

Figure 4. Temperature Deviation for Storm (3 hr Half-Life at z = 1 H) 
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